%\documentclass{cumcmthesis}
\documentclass[withoutpreface,bwprint]{cumcmthesis} %去掉封面与编号页，电子版提交的时候使用。


\usepackage[framemethod=TikZ]{mdframed}
\usepackage{url}   % 网页链接
\usepackage{subcaption} % 子标题
\usepackage{algorithm} % 算法
\usepackage{algpseudocode} % 伪代码
% 需在导言区加载必要宏包（若未加载）
\usepackage{booktabs}  % 用于绘制规范三线表
\usepackage{array}     % 用于灵活定义列格式
\usepackage{diagbox} % 绘制斜线表头
% \usepackage{tikz}
% \usetikzlibrary{shapes.geometric, arrows} % 用于图形和箭头

% % 定义流程图样式
% \tikzstyle{startstop} = [rectangle, rounded corners, minimum width=3cm, minimum height=1cm,text centered, draw=black, fill=red!30]
% \tikzstyle{process} = [rectangle, minimum width=3cm, minimum height=1cm, text centered, draw=black, fill=blue!30]
% \tikzstyle{decision} = [diamond, minimum width=3cm, minimum height=1cm, text centered, draw=black, fill=green!30]
% \tikzstyle{arrow} = [thick,->,>=stealth]

\title{基于目标优化构建烟幕干扰弹投放策略的研究}
\tihao{A}
\baominghao{xxxx}
\schoolname{XX大学}
\membera{ }
\memberb{ }
\memberc{ }
\supervisor{ }
\yearinput{2023}
\monthinput{09}
\dayinput{04}
\setlength{\parindent}{2em} % 设置首行缩进

\begin{document}

 \maketitle
 \begin{abstract}
    在现代作战体系中，烟幕干扰弹是一种有效对抗敌方制导武器的手段。在给定无人机和导弹初始位置及运动控制的条件下，本文针对烟幕干扰弹的投放策略展开研究，以最大化本方目标的遮蔽效果。

    \textbf{针对问题一}，我们根据题意所给投放策略，建立了导弹、无人机、烟幕弹及烟幕云团的运动学模型，给出上述物体随时间的位置坐标变化关系，
    并设计了基于数值模拟的求解流程，计算得出\textbf{导弹$M_1$被有效遮蔽的时长为1.39s}。

    \textbf{针对问题二}，我们构建最大化导弹$M_1$遮蔽时长的\textbf{单目标优化模型}，以无人机$FY_1$的飞行角度、速度、烟幕弹的投放时间和起爆时间作为决策变量，对飞行角度、速度等变量
    进行约束。由于决策变量数量较多，我们采用粒子群算法对问题进行求解，得到最优投放策略使得\textbf{导弹$M_1$被遮蔽时长达到4.60s}。

    \textbf{针对问题三}，我们在问题二的基础上增加了两枚烟幕弹的投放时间和起爆时间作为决策变量，建立了单目标优化模型。在尝试粒子群算法后，我们发现决策变量过多导致粒子群算法收敛性变差。
    因此，我们采用降维的思路，对各变量进行分步优化。首先得到最佳飞行角度$\theta=180^ \circ$，之后对飞行速度$v$进行变步长搜索，利用模拟退火算法对其余变量进行分步优化，
    最终得到三枚烟幕弹\textbf{对$M_1$的最优遮蔽时间为6.16s}。

    \textbf{针对问题四}，由于单目标优化模型的决策变量增加至12个，我们同样基于降维的思路，将问题四分拆为两个部分进行求解：首先，建立\textbf{航迹选择子模型}，根据导弹轨迹生成一系列起爆候选
    参考点，并将遮蔽效果转化为\textbf{球面投影覆盖}问题，从而简化对候选点的评估。根据评估得到的最佳候选点反解出最佳飞行角度$\theta^ *$及对应飞行速度$v^*$；其次，固定每架无人机的航向，构建降维后的单目标优化模型，
    求解过程中，我们在预估飞行速度附近对速度进行遍历搜索，并使用模拟退火算法对其余变量进行优化求解，最终得到三架无人机各投放一枚烟幕弹\textbf{对$M_1$的最优遮蔽时间为11.44s}。

    \textbf{针对问题五}，基于前四问的分析，我们首先通过枚举法求出各个无人机干扰各个导弹的各个最佳飞行方向（共15种情况），在此基础上\textbf{三级分层优化搜索}每架无人机的最佳速度，最终通过\textbf{模拟退火算法}优化单枚烟幕干扰弹的投放时间和引信延时，之后再进行\textbf{网格精细搜索}，在最优解附近的小范围内进行密集采样
    ， 得到\textbf{三枚导弹最优总体遮蔽时间为24.73s}。

\keywords{单目标优化模型\quad  粒子群算法\quad 网格搜索\quad 航迹选择}
\end{abstract}



\section{问题重述}
    \subsection{问题背景}
烟幕干扰弹通过在空气中施放大量的气溶胶微粒，以改变电磁波在介质中的传输特性，
从而可以有效地干扰和欺骗对方的光电系统和制导系统，既能隐藏真目标，
也能突显出假目标迷惑地方\cite{ref1}。随着烟幕干扰技术的发展，新一代的系统不仅
能将其精确抛撒至预定位置，还能精确控制起爆时间。在此背景下，投放策略的优劣将成为
影响遮蔽效果的关键因素。
\subsection{问题要求}
在该模型中，导弹速度为300m/s，直指在坐标原点的假目标，真目标为半径7m，高10m的圆柱体，其下底面
圆心位于(0,200,0)。雷达发现导弹时，3枚导弹$M_1$、$M_2$、$M_3$和5架无人机 $FY_1 - FY_5$的坐标均已知。无人机受领任务后瞬时设定航向和速度，
并立即匀速直线飞行。烟幕干扰弹脱离无人机后在重力作用下作平抛运动。在设定时间起爆后，烟幕弹形成遮蔽半径10m、有效遮蔽时长20s的烟幕云团。
并以3m/s的速度匀速下降，每架无人机投放相邻两枚烟幕干扰弹至少间隔1s。

\textbf{问题1}  
利用无人机$FY_1$投放1枚烟幕干扰弹实施对$M_1$的干扰，设定$FY_1$以120m/s的速度向假目标（原点）方向飞行，
受领任务1.5s后投放1枚烟幕干扰弹，间隔3.6s后起爆,据此计算烟幕干扰弹对M1的有效遮蔽时长。

\textbf{问题2}  
在问题一的基础上，本问不再给出$FY_1$的飞行方向、飞行速度、
烟幕干扰弹投放点位（时间）、烟幕干扰弹起爆点位（时间）等参数，需要对上述参数进行优化并给出最优解，使得干扰弹对$M_1$的遮蔽时间最大化。

\textbf{问题3} 
在问题二的基础上，无人机$FY_1$投放烟幕干扰弹的数量增加为3枚，要求确定$FY_1$的飞行方向、飞行速度和3枚烟幕干扰弹的各自
投放点位（时间）、起爆点位（时间），使得$M_1$的总体被遮蔽时间最长。

\textbf{问题4}  
利用$FY_1$、$FY_2$、$FY_3$这3架无人机各投放1枚烟幕干扰弹来对导弹$M_1$进行遮蔽，给出这3架无人机各自的飞行方向和速度，以及
3枚烟幕干扰弹各自的投放点位（时间）、起爆点位（时间），使得$M_1$的总体被遮蔽时间最长。

\textbf{问题5}
利用全部5架无人机，且每架无人机至多投放3枚烟幕干扰弹，来对全部3枚导弹进行遮蔽，要求给出对五架无人机分别给出最优投放
策略（飞行方向和速度，投放和起爆点位），使得三枚导弹的总体被遮蔽效果最佳。

\section{问题分析}
    首先明确达到遮蔽效果的条件：若烟幕云团在有效时间内与导弹位置到真目标圆柱体任意一点的
    有限线段相交，则视为达到遮蔽效果。以此为基础，对问题展开以下分析：
    \subsection{问题一分析}
     问题一要求我们计算在给定条件下烟幕干扰弹对$M_1$的有效遮蔽时长。首先，我们对导弹、无人机、烟幕弹
     以及云团进行运动学模型刻画；接着，将起爆后的20s的时间段以0.01s为步长遍历，根据运动学模型求解每个时刻各个物体的位置坐标，
     并判定该时刻云团是否达到遮蔽效果；最后统计所有被遮蔽的时刻，即可得到所求。
    \subsection{问题二分析}
    问题二要求我们设计无人机$FY_1$的投放策略，使得导弹$M_1$在20s内被有效遮蔽的时长最大。我们可以将问题二建模为一个单目标优化模型，
    以无人机的飞行角度、速度、烟幕弹的投放时间和起爆时间作为决策变量，导弹$M_1$的遮蔽时长作为目标函数，并根据题意构建约束条件。
    由于决策变量较多且目标函数为非光滑分段函数，我们采用粒子群优化算法（PSO）对问题二进行求解。
    \subsection{问题三分析}
    与问题二相比，问题三无人机$FY_1$投放烟幕干扰弹的数量增加为3枚。我们对问题二的模型进行扩展，增加两枚烟幕弹的投放时间和起爆时间作为决策变量，
    并构建类似的单目标优化模型。由于决策变量进一步增多，我们在尝试粒子群算法的基础上，结合模拟退火算法对各个决策变量进行分步优化，从而降低每一步优化求解的维度。
    \subsection{问题四分析}
    问题四将无人机数量增加到3架，决策变量的数量增至12个，进一步增加了模型的复杂度。为更加有效地降维，我们对问题四分为以下两步求解：首先根据几何投影评估建立
    航际选择模型求解每架无人机的最佳飞行方向，并给出预估的飞行速度，在此基础上，结合问题三的思路，在预估飞行速度附近进行遍历搜索，并使用模拟退火算法对其余变量进行优化求解。
    \subsection{问题五分析}
    基于前四问的分析，我们在问题五采用分布优化思路求解是首先通过枚举法求出各个无人机干扰各个导弹的各个最佳飞行方向（共15种情况），在此基础上遍历搜索每架无人机的最佳速度，最终通过模拟退火算法优化单枚烟幕干扰弹的投放时间和引信延时，之后再进行网格精细搜索，在最优解附近的小范围内进行密集采样。


\section{模型假设}
    \begin{enumerate}
        \item 除真实目标外，忽略所有无人机、导弹等物体的体积，视为质点处理。
        \item 忽略烟幕浓度衰减对遮蔽效果的影响，将有效遮蔽范围简化为固定半径10m，以云团中心为圆心的球形区域。
        \item 忽略气象条件如风速、风向、能见度等对无人机飞行、烟幕弹运动和烟幕云团遮挡效果的影响。
        \item 假设在飞行过程中，任何无人机、导弹间不会发生相对碰撞。
        \item 假设导弹雷达的探测范围为360$^\circ$的无限远空间，除烟幕云团外不存在其他任何遮挡。
    \end{enumerate}

\section{符号说明}
    \begin{table}[!htbp]
    \caption{符号说明}\label{tab:001} \centering
    \begin{tabular}{ccccc}
        \toprule[1.5pt]
        符号 & 单位 & 说明 \\
        \midrule[1pt]
        $\theta$ & $^\circ$ & 无人机飞行角度，与x轴正方向夹角 \\
        $v$ & m/s & 无人机飞行速度 \\
        $td$ & s & 烟幕弹投放时间 \\
        $te$ & s & 烟幕弹起爆时间 \\
        $d_{PM \rightarrow C}$ & m & 真目标-导弹连线$PM$到烟幕球中心$C$的距离 \\
        $O_k(t)$ & 0-1指示量 & t时刻对导弹$M_k$的遮蔽状态 \\
        $T$ & s & 有效遮蔽总时长 \\
        $S(t)$ & - & 瞬时遮蔽率 \\
        \bottomrule[1.5pt]
    \end{tabular}
\end{table}

\section{模型准备}

    \subsection{对导弹运动的刻画} 
    根据题意，三枚导弹从各自初始位置指向假目标（即原点）以恒定速度$v_m=300m/s$飞行，导弹的位置为关于时刻的函数。
    设导弹$M_k$在时刻$t$的位置为$m_k(t)$，由题干条件，三枚导弹的运动方向单位向量为：
    $$
        \vec{n_1}=\left(-0.9950,0,-0.0995\right), \quad
        \vec{n_2}=\left(-0.9934,-0.0314,-0.1098\right), \quad
        \vec{n_3}=\left(-0.9939,0.0331,-0.1049\right). \quad
    $$
    由此可得导弹的位置随时间变化关系如下:
    $$
        m_k(t)=m_k(0)+v_m\cdot t\cdot \vec{n_k}, \quad k=1,2,3.
    $$
     
    \subsection{对无人机运动的刻画} 
    根据题意，无人机$FY_i$从各自初始位置以指定速度$v_i \left( v_i \in [70, 140] \right)$沿指定方向飞行，设无人机$FY_i$在时刻$t$的位置为$u_i(t)$。
    题干限制无人机的飞行方向均在水平面内，设无人机$FY_i$的飞行方向与x轴正方向夹角为$\theta_i$，则无人机的位置随时间变化函数可表示如下：
    $$
        u_i(t)=u_i(0)+v_i\cdot t\cdot \left(\cos\theta_i,\sin\theta_i,0\right), \; i=1,2,3,4,5.
    $$
    
    \subsection{对烟幕弹及烟幕云团运动的刻画} 
    记无人机$FY_i$投放的第j颗烟幕弹为$C_{ij}$，其在时刻$t$的位置为$c_{ij}(t)$，投放时间为$td_{ij}$，起爆时间为$te_{ij}$。

    烟幕弹在投放瞬间位置坐标即为无人机坐标$u_{ij}(t)$，在投放后$td_{ij}$到$te_{ij}$时间段内，水平方向速度保持为$v_i$不变，竖直方向受重力加速度作用，做平抛运动。起爆时间$te_{ij}$瞬间扩散
    为烟幕球，并以$v_{sink} = 3m/s$的速度下沉，烟幕球的半径为$r=10m$，有效时间区间为$[te_{ij}, te_{ij}+20]$。烟幕弹（烟幕球中心）位置坐标随时间变化关系如下:
    $$
        s_{ij}(t)=
        \begin{cases}
            u_i(td_{ij})+v_i\cdot (t-td_{ij})\cdot (\cos\theta_i,\sin\theta_i,0)-\left(0,0,\frac{1}{2}g(t-td_{ij})^2\right), & t \in [td_{ij},te_{ij}),\\
            c_{ij}(te_{ij})-\left(0,0,v_{sink}\cdot (t-te_{ij})\right), & t \in [te_{ij},te_{ij}+20].
        \end{cases}
    $$
    \subsection{对有效遮蔽范围的判定}
    将导弹的侦测过程抽象为与真目标的直线连线，烟幕云团在起爆有效时间范围$[te_{ij},te_{ij}+20]$内的遮蔽范围为半径10m的球形区域。

    \begin{figure}
    \centering
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=.95\textwidth]{A1_initial.png}
        \caption{初始状态各元素相对位置图}
        \label{fig:initial}
    \end{minipage}
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=.95\textwidth]{pic5.4.png}
        \caption{烟幕有效遮蔽导弹示意图}
        \label{图5.4}
    \end{minipage}
    \caption{模型空间关系示意图}
    \label{fig:pro_prepare}
    \end{figure}

    则“有效遮蔽”的判定可表述如下：

    对于导弹$M$到真目标上任意一点$P$的连线$PM$，必须在t时刻与烟幕球存在交点，则导弹$M$在t时刻是被有效遮蔽的。

    真目标为底面中心$H$位于坐标(0,200,0)，半径$r=7m$,高$h=10m$的圆柱体。由于有效遮蔽的判定条件需要求连线$PM$中距离烟幕球球心最大值，
    我们考虑在圆柱顶面圆周和底面圆周的各点进行均匀采样n个点，作为真目标上点$P$的候选集合，如\cref{图5.4}所示。

    设采样点集合为$\{P_i\}, i=1,2,\dots,2n$，其中前n个点为底面圆周采样点，后n个点为顶面圆周采样点，则各点$P_i$的坐标为：
    $$
        P_i=
        \begin{cases}
            (r\cos\frac{2\pi i}{n}, 200 + r\sin\frac{2\pi i}{n}, 0), & i=1,2,\dots,n,\\
            (r\cos\frac{2\pi (i-n)}{n}, 200+r\sin\frac{2\pi (i-n)}{n}, h), & i=n+1,n+2,\dots,2n.
        \end{cases}
    $$

    设t时刻烟幕云团中心$C$与连线$PM$的距离为$d_{PM \rightarrow C}(t)$。注意此处的距离并不能简单判定为点到直线的距离，
    而是云团中心到以$P, M$为端点的线段的距离最小值，我们使用投影法进行处理：
    $$
        \overrightarrow{MP}=(P_x, P_y, P_z)-(M_x, M_y, M_z), \quad \overrightarrow{MC}=(C_x, C_y, C_z)-(M_x, M_y, M_z).
    $$
    则线段$MP$上点的参数式可表示为$X(u) = M + u \cdot \overrightarrow{MP}, \; u \in [0,1]$。将$C$正交投影到直线$PM$上，则有：
    $$
        u_C=\frac{\overrightarrow{MP} \cdot \overrightarrow{MC}}{\left| \overrightarrow{MP} \right|^2}.
    $$
    $u_C$的取值范围决定投影点$Q = M + u_C \cdot \overrightarrow{MP}$的位置，由此我们可以计算得到点$C$到线段$PM$的距离：
    $$
        d_{PM \rightarrow C}(t)=
        \begin{cases}
        \left| \overrightarrow{MP} \right|, & u_C \leq 0, \text{投影点在} M \text{外侧}\\
        \left| \overrightarrow{CQ} \right|, & 0 < u_C < 1, \text{投影点在} PM \text{线段上}\\
        \left| \overrightarrow{MC} \right|, & u_C \ge 1, \text{投影点在} P \text{外侧}
        \end{cases}
    $$
    遍历整个真目标上的任意一点P，当$\max{d_{PM \rightarrow C}(t)} \leq 10$时，即为有效遮蔽状态。记遮蔽指示量$O_k(t)$为t时刻对导弹$M_k$的遮蔽状态，若被遮蔽则$O_k(t)=1$，否则$O_k(t)=0$。$O_k(t)$的判定方式如下：\\
    $$
        O_k(t)=
        \begin{cases}
            1, & \exists C_{ij}\; s.t. \; \max{d_{PM \rightarrow C}(t)}  \leq 10 \; t\in[te_{ij},te_{ij}+20],\\
            0, & Otherwise.
        \end{cases}
    $$

\section{模型的建立与求解}
    \subsection{问题一模型的求解}
        问题一给出了一个具体的烟幕弹投放策略，包括无人机$FY_1$的具体飞行方向、速度、投放时间和起爆时间，要求干扰弹对$M_1$的有效遮蔽时长。
        可通过仿真模拟，直接根据上述已经建立好的基本模型对其进行求解。

        通过题目所给数据，我们易求得烟幕弹投放和起爆时的位置坐标及起爆后任意时刻烟幕球中心坐标$c_{11}(t)$，导弹位置$m_1(t)$已知。以起爆时间为开始时间，将
        [0,20]这个时间段以0.01s为步长遍历。对于每个时刻，我们对真目标上下边缘各均匀抽取4个点位，计算每个点位与导弹$M_1$的连线与烟幕球中心的距离，取其中最大值
        判定该时刻导弹$M_1$是否被遮蔽。最终统计所有时刻中被遮蔽的时长即为所求。
        \begin{figure}[!h]
        \centering
        \includegraphics[width=0.8\textwidth]{model1_flow.png}
        \caption{问题一求解流程图}
        \label{fig:pro1_flow}
        \end{figure}

        通过编程实现上述流程（代码见附录B），计算得出该投放策略下，导弹$M_1$被有效遮蔽的时长为1.39s，有效遮蔽时间区间为起爆时刻后[2.96, 4.35](s)。
        取时间间隔为3秒，0到20秒，各时间点导弹$M_1$与真目标连线到烟幕球中心的距离变化情况如\cref{fig:pro1_result}所示：
        
    

    \begin{figure}
    \centering
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=0.95\textwidth]{problem1_1.png}
        \subcaption{[0,20]s全局图}
        \label{fig:pro1_result_a}
    \end{minipage}
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=0.95\textwidth]{problem1_2.png}
        \subcaption{有效遮蔽区间局部放大图}
        \label{fig:pro1_result_b}
    \end{minipage}
    \caption{问题一烟幕球中心距离变化图}
    \label{fig:pro1_result}
    \end{figure}

    \textbf{结果分析：}
    
    由\cref{fig:pro1_result_a}可见，在起爆后[0,5](s)内，烟幕球中心与连线距离基本保持不变，维持在一个较小的范围内，而在5s之后距离以较快的速度增大，遮蔽效果迅速消失。

    由局部放大图\cref{fig:pro1_result_b}可见，在起爆后4.4s左右，距离函数在平缓降至最小值后发生突变，迅速增大。为分析该突变产生原因，我们提取问题一求解过程中各时间点的
    导弹坐标$m_1(t)$、烟幕球中心坐标$c_{11}(t)$，以渐变色表示时间先后，绘制导弹和烟幕球轨迹如\cref{fig:pro1_visual_a}。由问题一条件限制，导弹和无人机运动轨迹均在$xOz$平面内，
    故\cref{fig:pro1_visual_a}将三维坐标投影到$xOz$平面后简化为二维轨迹处理。

    作几何示意图如\cref{fig:pro1_visual_b}，在$t=4.4s$时刻前，烟幕球中心在直线上的投影点$Q$位于线段$PM$上，距离函数$d_{PM \rightarrow C}(t) = \left| \overrightarrow{CQ} \right|$，
    其主要变化因素为烟幕球中心$C$的下沉运动，$d_{PM \rightarrow C}(t) = \left| \overrightarrow{CQ} \right|$变化相对平稳。而在$t=4.4s$后，投影点$Q$移动到线段$PM$的端点$P$外侧，
    $d_{PM \rightarrow C}(t) = \left| \overrightarrow{CQ} \right|$表达式变为$\left| \overrightarrow{MC} \right|$，且受导弹移动影响为主。由于导弹移动速度$v_0=300m/s$远大于云团
    下降速度，$M$点位置迅速变化，使$d_{PM \rightarrow C}(t) = \left| \overrightarrow{CQ} \right|$在最小值附近出现明显拐点。

    \begin{figure}
    \centering
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=0.95\textwidth]{problem1_trajectory.png}
        \subcaption{问题一轨迹可视化}
        \label{fig:pro1_visual_a}
    \end{minipage}
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=0.95\textwidth]{pro1_geo.pdf}
        \subcaption{距离函数突变几何示意图}
        \label{fig:pro1_visual_b}
    \end{minipage}
    \caption{问题一距离突变分析示意图}
    \label{fig:pro1_visual}
    \end{figure}
 \textbf{若将真目标简化为圆柱体底面中心的一个质点}，采用相同的方法可以计算得出有效遮蔽区间[2.94, 4.35](s)，总有效遮蔽时长约为1.41s。
    将真目标视为质点与多点采样的结果非常接近，相对误差仅为约1.4\%。这表明在当前参数下，是否考虑真目标本身的形状大小对结果影响较小，在问题三及之后模型的建立和求解过程当中，我们考虑在优化投放参数时忽略真目标形状，将其视为圆柱体底面中心的质点以简化计算，
    而在实际求解有效遮蔽时间是仍将目标视为圆柱体处理，其合理性分析如下：

    \begin{enumerate}
        \item 从几何因素考虑，烟幕云团是一个半径为10m的球体，而圆柱体目标的尺寸相对较小（半径7m、高10m）。云团足够大，能够同时覆盖圆柱体的中心点和边缘点的大部分区域。因此，检查中心点是否被遮蔽可以近似代表整个目标被遮蔽的情况，误差在可接受范围内。
        \item 从计算效率考虑，在后续问题3、4、5中，涉及多个烟幕弹和多个导弹，计算复杂度会显著增加。使用质点模型可以大大简化计算，提高优化算法的效率，而结果仍能保持较高的准确性。
        \item 在精度和计算成本之间权衡，既然质点模型与代表点模型的结果接近，说明质点模型足以捕捉主要物理现象，这样简化有利于快速迭代和求解多变量优化问题。
    \end{enumerate}

    %\newpage
    \subsection{问题二模型的建立与求解}
    \subsubsection{模型建立}
    问题二仍然为无人机$FY_1$投放烟幕弹对导弹$M_1$的遮蔽优化问题，并要求我们设计具体的投放策略，即给出无人机的飞行方向$\theta$、速度$v$、投放时间$td$和起爆时间$te$，使得导弹$M_1$在20s内被有效遮蔽的时长最大。

    由模型准备部分，给定时间t，导弹$M_1$位置$m_1(t)$已知，烟幕球中心位置$c_{11}(t)$可由$\theta, v, td, te$计算得到，进而可计算出$t$时刻导弹$M_1$是否被遮蔽的指示量$O(t;\theta,v,td,te)$。

    定义决策向量
    $$
    x = (\theta, v, t_d, t_e)
    $$
    目标函数为遮蔽时间函数$T(x)$
    $$
    T(x) = \int_{te}^{te+20} O(t; x) dt
    $$
    为便于程序仿真，对时间以0.1s为步长进行离散处理，通过近似求和代替积分，遮蔽时间函数可表示为：
    $$
    T(x) = \sum_n O(t_n; x) \Delta t \quad t_n \in [te, te+20], \Delta t = 0.1s
    $$
    其中遮蔽状态函数$O_k(t_n)$沿用模型准备中的几何判定方式，即：
    $$
    O_k(t_n)=
            \begin{cases}
                1, & \exists C_{ij}\; s.t. \; \max{d_{PM \rightarrow C}(t_n)}  \leq 10 ,\\
                0, & Otherwise.
            \end{cases} \\
    $$
    最终得到的模型如下所示：
    \begin{align*}
            \max\limits_{(\theta, v, t_d, t_e)} T(x) &= \sum_n O(t_n; x) \Delta t \quad t_n \in [te, te+20]\\
            s.t. \begin{cases}
                & \theta  \in [0, 2\pi], \\
                & v \in [70, 140], \\
                & td \leq te. \\
            \end{cases}
    \end{align*}

    \subsubsection{模型求解}
    由于本问题中目标函数为非光滑分段函数，难以通过解析方式求解。我们采用该算法PSO粒子群优化算法对问题二进行求解。该算法具有收敛速度快、全局搜索能力强等优点，其主要实现思路如\cref{alg:PSO_pro2}.
 

\begin{algorithm}
\caption{问题二PSO优化算法}
\label{alg:PSO_pro2}
\begin{algorithmic}[1]
\Require 粒子数 $N$，最大迭代 $K$，惯性权重 $\omega$，认知系数 $c_1$，社会系数 $c_2$，边界 $\mathbf{l},\mathbf{u}$，命中时间 $t_{\mathrm{hit}}$
\Ensure 近似全局最优解 $\mathbf{x}_{\mathrm{PSO}}$
\State 随机初始化粒子位置 $\mathbf{x}_i^{(0)}\sim \mathcal{U}(\mathbf{l},\mathbf{u})$，速度 $\mathbf{v}_i^{(0)}$，$i=1,\dots,N$
\State 对每个粒子做一次可行性修复
\State 计算 $\tilde f(\mathbf{x}_i^{(0)})$，设置个体最优 $\mathbf{p}_i^{(0)}\leftarrow\mathbf{x}_i^{(0)}$，群体最优 $\mathbf{g}^{(0)}\leftarrow\arg\min_i \tilde f(\mathbf{x}_i^{(0)})$
\For{$k=0$ \textbf{to} $K-1$}
    \For{$i=1$ \textbf{to} $N$}
        \State 采样 $r_1,r_2\sim\mathcal{U}(0,1)$
        \State 更新速度$\mathbf{v}_i^{(k+1)}=\omega\,\mathbf{v}_i^{(k)}
        +c_1 r_1(\mathbf{p}_i^{(k)}-\mathbf{x}_i^{(k)})
        +c_2 r_2(\mathbf{g}^{(k)}-\mathbf{x}_i^{(k)})
        $
        \State 更新位置 $\mathbf{x}_i^{(k+1)}=\mathbf{x}_i^{(k)}+\mathbf{v}_i^{(k+1)}$
        \State 边界投影 $\mathbf{x}_i^{(k+1)}\leftarrow \min\{\mathbf{u},\,\max(\mathbf{l},\,\mathbf{x}_i^{(k+1)})\}$
        \State 可行性修复
        \State 计算适应度 $\tilde f(\mathbf{x}_i^{(k+1)})$
        \If{$\tilde f(\mathbf{x}_i^{(k+1)})<\tilde f(\mathbf{p}_i^{(k)})$} \State $\mathbf{p}_i^{(k+1)}\leftarrow\mathbf{x}_i^{(k+1)}$
        \Else \State $\mathbf{p}_i^{(k+1)}\leftarrow\mathbf{p}_i^{(k)}$
        \EndIf
    \EndFor
    \State 更新群体最优 $\mathbf{g}^{(k+1)}\leftarrow\arg\min_{i}\tilde f(\mathbf{p}_i^{(k+1)})$
    \State （可选）调整 $\omega$ 或 $c_1,c_2$
\EndFor
\State \Return $\mathbf{x}_{\mathrm{PSO}}=\mathbf{g}^{(K)}$
\end{algorithmic}
\end{algorithm}

    为保证算法效率，在PSO算法中适当增大时间步长$\Delta t=0.5s$，得到粗糙的最优解$\mathbf{x}_{\mathrm{PSO}}$后，再以该解为初始点，对局部空间进行搜索，得到最终的全局最优解$x_{\mathrm{final}}$，并缩减时间步长$\Delta t$ 为0.1s，得到精确的有效遮蔽时长。


    取粒子数$N=50$，最大迭代次数$K=100$，认知系数和社会系数分别为$c_1=c_2=2.0$，边界$\mathbf{l}=(0,70,0,0), \mathbf{u}=(2\pi,140,20,20)$，并通过编程实现上述算法（代码见附录B）.求解结果如下：
    $$
    x_{\mathrm{PSO}} = (\theta, v, t_d, t_e) = (7.18^\circ, 99.99m/s, 0.0s, 0.59s)
    $$
    对应的投放点和起爆点坐标分别为
    \begin{align*}
        (X_d, Y_d, Z_d) = (17800, 0, 1800), \\
        (X_e, Y_e, Z_e) = (17858.19, 7.33, 1798.31)
    \end{align*}

    该投放策略下，导弹$M_1$被有效遮蔽的时长为4.60s，有效遮蔽时间区间为起爆时刻后[2.55, 7.15](s)，云团中心到导弹-真目标连线距离随时间变化情况如\cref{fig:pro2_result}所示。

     
    \begin{figure}
    \centering
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=0.95\textwidth]{problem2_1.png}
        \subcaption{[0,20]s全局图}
        \label{fig:pro2_result_a}
    \end{minipage}
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=0.95\textwidth]{problem2_2.png}
        \subcaption{有效遮蔽区间局部放大图}
        \label{fig:pro2_result_b}
    \end{minipage}
    \caption{问题二烟幕球中心距离变化图}
    \label{fig:pro2_result}
    \end{figure}

   

    \subsection{问题三模型的建立与求解}
    \subsubsection{模型建立}
    问题三在问题二的基础上将无人机$FY_1$投放烟幕弹的数量增加到3枚，同样要求设计具体的投放策略。根据问题一、二的分析，我们建立如下单目标优化模型：
    
    \textbf{目标函数：}
    $$
        \max\limits_{x} T(x) = \sum_n O(t_n; x) , \Delta t = 0.1s
    $$

    \textbf{决策变量$x$与约束条件：}
    $$
        x = (\theta, v, td_{i}, te_{i})
    $$
    其中，$\theta, v$为无人机飞行方向和速度，$td_{i}, te_{i}$分别为第i枚烟幕弹的投放时间和起爆时间，$i=1,2,3$。
    由题干限制，无人机投放相邻两枚烟幕弹的时间间隔至少为1s，故由以下约束条件：
    \begin{align*}
        td_{2}-td_{1} \ge 1, \\
        td_{3}-td_{2} \ge 1.
    \end{align*}
    根据模型准备中对导弹位置$m_1(t)$和烟幕球位置$c_{ij}(t)$的刻画，以及遮蔽状态函数$O(t_n; x)$的几何判定方式，最终得到的模型如下所示：
    \begin{align*}
            \max\limits_{x=(\theta, v, t_{di}, t_{ei})} T(x) = \sum_n O(t_n; x) \Delta t \\
            s.t. \begin{cases}
               & \theta \in [0, 2\pi], \\
               & v \in [70, 140], \\
               & td_{1} \leq td_{2}-  1 \leq td_{3}-2, \\
               & td_{i} \leq te_{i}, i=1,2,3. \\
            \end{cases}
    \end{align*}

    \subsubsection{模型求解}
    尝试采用分布式粒子群算法进行求解后，我们发现，由于变量维度较高，粒子群算法在该问题上表现出较差的收敛性，难以得到满意的结果。
    因此我们考虑对问题三进行降维处理，将整体优化转变为对角度、速度和投放/起爆时间的分步优化。

    通过修改粒子群优化算法并对$0^\circ$到$360^\circ$进行搜索，我们发现，在$180^\circ$附近
    才能让3枚烟幕弹产生的遮蔽效果较为均衡，而其他角度下最多只有一枚烟幕弹起作用。
    考虑固定$\theta = 180^\circ$，从而问题三的决策变量为：
    $$
        x = (v, t_{di}, t_{ei})
    $$
    为进一步降维，我们对飞行速度$v$在[70,140]区间进行变步长搜索找出最优速度值，从而求出整体优化策略。取均匀分布的n个速度值，对每个$v_j(j=1,2,...,n)$，通过模拟退火算法对投放时间$t_{di}$和起爆时间$t_{ei}$进行联合优化，得到
    每个飞行速度$v_j$对应的最优投放策略及其有效遮蔽时长$T_j$，并通过序列化优化方法依次对3枚烟幕弹进行优化，提升了算法效率。

    根据代码求解结果（代码见附录B），速度在121m/s附近所得的总遮蔽时间呈现出最优值，如\cref{fig:shielding_time3}所示。
    

    \begin{figure}[!h]
    \centering
     \includegraphics[width=.8\textwidth]{shielding_time.png}
     \caption{遮蔽时间对比图}
     \label{fig:shielding_time3}
 \end{figure}
    所以，我们继续对121m/s的邻域取更小的步长($\Delta t = 0.02s$)来进一步搜索。最终，我们得到，在速度
    v=121.10m/s时，总遮蔽时长最大，为6.16s。各烟幕弹遮蔽时间对比如\cref{tab:small}所示。
 
    \begin{table}[H]
        \caption[标签名]{小步长遮蔽时间对比}\centering
        \label{tab:small}
        \begin{tabular}{ccccc}
            \toprule[1.5pt]
            速度v(m/s) & 总遮蔽时间/s&烟幕弹1遮蔽时间/s  & 烟幕弹2遮蔽时间/s&烟幕弹3遮蔽时间/s  \\
            \midrule[1pt]
            121.00 & 6.12 &3.08 &2.14 &0.94\\
            121.02 & 6.12 &3.08 &2.14 &0.94\\
            121.04 & 6.12 &3.08 &2.14 &0.94\\
            121.06 & 6.12 &3.08 &2.14 &0.94\\
            121.08 &6.12  &3.08 &2.14 &0.94\\
            121.10 &6.16 &3.08  &2.20 &0.94\\
            121.12 &5.16 &3.02 &1.56 & 0.76\\
            121.14 & 5.14 &3.02 &1.62 &0.76\\
            
            \bottomrule[1.5pt]
        \end{tabular}
    \end{table}

    同时，我们得到了飞行角度为180$^\circ$，速度为121.10m/s时3枚烟雾弹的投放、引爆时间以及遮蔽时间的区间。如\cref{tab:data3}所示。

    \begin{table}[H]
        \caption[short]{各烟幕弹数据}\centering
        \label{tab:data3}
        \begin{tabular}{ccccc}
            \toprule[1.5pt]
            烟幕弹序号& 投放点坐标 & 起爆点坐标 & 有效遮蔽时长(s) &有效遮蔽区间\\
            \midrule[1pt]
            1 & (17800, 0,1800) &(17401, 0, 1747.03) & 3.08 &[5.64, 8.72]\\
            2 & (17362, 0, 1800) &(16742, 0, 1671.48) & 2.14 &[8.74, 10.94]\\
            3 & (17153, 0, 1800)&(16482.89, 0, 1649.84) & 0.94 &[10.88, 11.82]\\
            \bottomrule[1.5pt]
        \end{tabular}
    \end{table}

    \subsection{问题四模型的建立与求解}
    问题四中无人机数量增加到3架，每架无人机的投放烟幕弹数量限制在1枚，针对导弹$M_1$进行有效遮蔽时长最大化。如果同前述问题一样建立单目标优化模型，则决策变量将增加到12个，求解难度和复杂度大大增加。因此我们借鉴问题三的降维思路，将问题四的解决分为两步：首先，
    建立如下的航迹选择模型，以较小的时间代价求得每个无人机$FY_i$的最佳航向$\theta^*_i$，并给出参考起爆位置和无人机飞行速度$v_i^*$；之后，固定每架无人机的航向，调用问题三的模型优化其余变量，得到最佳调度方案。

    \subsubsection{航迹选择子模型建立}
    由于前述问题的模型在评估投放策略时需数值模拟云团和导弹的运动过程，求解较为复杂。本模型提出一种较为简便的几何投影方法用于评估。
    
    在任意时刻$t$，导弹位置$m(t)$已知，导弹视觉范围可表示为以$m(t)$为球心的单位球$\mathcal{S}^2$内，则真目标和云团在该视觉球内投影的重合部分占比越大，说明云团对导弹遮蔽效果越好。
    
    记真目标上所有点集合为$\mathcal{P}$，假设已知该时刻云团位置$C$.从导弹出发的视线可用单位向量$u \in \mathcal{S}^2$表示，即射线：
    $$
        L(u) = \{ m(t) + s u | s \ge 0 \}
    $$
    则真目标在导弹的视觉球内的投影域为：
    $$
    \mathcal{D}_P(t) = \{ u \in \mathcal{S}^2 | L(u) \cap \mathcal{P} \neq \emptyset \}.
    $$
    真目标为圆柱体，为方便计算，将$\mathcal{D}_P(t)$ 近似为视觉球上一球面圆，其中圆心为导弹-真目标中心连线$HM$与视觉球的交点，其角半径$\beta$可由下式计算得到：
    \begin{align*}
     \beta &= \arctan \left(\frac{R'}{|HM|}\right) \\
     R' &= \sqrt{R^2 + \left(\frac{h}{2}\right)^2}
    \end{align*}
    其中$R'$为圆柱体外接球的半径，$|HM|$为导弹与真目标中心的距离。

    同理，云团在导弹视觉球内的投影域也可用一球面圆$\mathcal{D}_C(t)$表示，记其角半径为$\alpha$：
    \begin{align*}
        \mathcal{D}_C(t) &= \{ u \in \mathcal{S}^2 |\angle (u, \overrightarrow{CM}) \leq \alpha \}, \\
        \alpha &= \arcsin \left(\frac{r}{|CM|}\right) 
    \end{align*}
    由此定义瞬时遮蔽率$S(t)$为：
    \begin{align*}
        S(t) &= \frac{Area(\mathcal{D}_P(t) \cap \mathcal{D}_C(t))}{Area(\mathcal{D}_P(t))} 
    \end{align*}
    考虑到真目标距离导弹较远，$\beta$数值较小，可以将球面圆重合问题近似转化为平面圆重合问题进行计算.记两圆圆心在视觉球上夹角为$\gamma$，则球面视觉域可对应为平面上两圆心相距$\gamma$，半径分别为$\alpha, \beta$的圆，遮蔽率计算如下：
    \begin{align*}
        S(t)&= \frac{\pi \min(\alpha^2, \beta^2) - \pi \left( \frac{\alpha^2 + \beta^2 - \gamma^2}{2} \right)}{\pi \beta^2} \\
        &= 1 - \frac{\alpha^2 + \beta^2 - \gamma^2}{2\beta^2}
    \end{align*}
    根据上述分析，通过计算遮蔽率$S(t)$可以较为简便地评估任意时刻$t$，为后续航向选择模型作准备。

    由题干条件，导弹命中假目标的时间$t_{hit} = 66.67s$，即云团能形成有效遮蔽的时间范围为$[0,t_{hit}]$。
    在$[0,t_{hit}]$上以0.5s为步长均匀采样得到时间点集$T_s= t_i, i=1,2,...,133$，可以计算出导弹在一系列时刻的位置坐标$m(t_i)$。对任意时刻$t_i$，
    导弹与真目标中心的连线$HM$已知，可通过$HM$确定云团起爆位置候选集$\mathcal{E} = E_1, E_2, ..., E_n$.候选集可通过一系列起爆高度平面$Z_E = z_1, z_2, ..., z_n$与$HM$交点得到：
    $$
        E_j = m(t_i) + \frac{z_j - z(m)}{z(H) - z(m)} \overrightarrow{HM}, \; j=1,2,...,n.
    $$
    其中n取120。为保证三架无人机的遮蔽时间区间尽量错开，我们限定每架无人机起爆高度在不低于其飞行高度以下1200m，即$z_j \in [z(U)-1200, z(U)]$，$z(U)$为无人机初始位置高度。
    
    建立评估函数$\varPi$如下：
    $$
        \varPi = S(t_i) + \lambda_v \frac{|E_{xy} - U_{xy}|}{t_i} - \lambda_t t_i - \lambda_{arg} \gamma + \lambda_z (z(U) - z(E))
    $$
    其中$S(t_i)$为该时刻的遮蔽率，$E_{xy}, U_{xy}$分别为起爆点和无人机初始位置在$xOy$平面的投影坐标，$z_U, z_E$分别为无人机初始位置和起爆点的高度，$\gamma$为云团起爆位置与导弹-真目标中心连线夹角，
    $\lambda_v, \lambda_t, \lambda_{arg}, \lambda_z$为权重系数。

    评估内容包括遮蔽率、预计飞行速度、飞行时间、航向偏离程度和高度差五个因素，并以遮蔽率$S(t_i)$为主导。经过多次调试，取$\lambda_v = 0.05, \lambda_t = 0.1, \lambda_{arg} = 0.2, \lambda_z = 0.15$，可较好反映
    起爆位置的遮蔽效果。

    遍历时间集$T_s$，对每一时间点$t_i$，计算其对应的起爆位置候选集$\mathcal{E}$，并计算评估函数$\varPi$。对每架无人机，取其在整个时间集$T_s$中$\varPi$降序排序前5位对应的起爆位置$E^*$，可分别反推出其对应的推荐航向$\theta^*$和飞行速度$v^*$为：
    \begin{align*}
        \theta^*_k &= \arctan \frac{y(E^*_k) - y(U)} {x(E^*_k) - x(U)}, \; k=1,2,3,4,5. \\
        v^*_k &= \frac{|E_{xy} - U_{xy}|}{t_i}
    \end{align*}
    \subsubsection{单目标优化子模型}
    由航向选择模型得到每个给定推荐航向$\theta_{ik}^*, i=1,2, \; k=1,2,3,4,5$.由于航向选择模型中对飞行速度、起爆时间等的计算较为粗糙，我们依旧建立
    单目标优化子模型对其余变量进行优化：
    \textbf{目标函数：}
    $$
        \max\limits_{x} T(x) = \sum_n O(t_n; x) , \Delta t = 0.1s
    $$
    对任意时刻$t$，函数$O(t; x)$需同时判断三个无人机投放的烟幕弹对导弹$M_1$的遮蔽状态，若任意一枚烟幕弹在$t$时刻对导弹$M_1$有效遮蔽，则$O(t; x)=1$，否则为0，即：
    $$
        O(t; x)=
            \begin{cases}
                1, & \exists C_{ii}\; s.t. \; \max{d_{PM \rightarrow C}(t)}  \leq 10 , \;i=1,2,3,\\
                0, & Otherwise.
            \end{cases} \\
    $$

    \textbf{决策变量$x$与约束条件：}
    航向选择模型对于每架无人机给出了5个推荐航向，我们固定每架无人机的航向$\theta_{ik}^*$，则问题四的决策变量为：
    \begin{align*}
        & x = (x_1, x_2, x_3), \\
        & x_i = (v_i, t_{di}, t_{ei}) \; i=1,2,3.
    \end{align*}
    无人机$FY_i$对应的决策变量为$x_i$，包括其飞行速度、投放时间和起爆时间。其中推荐飞行速度$v_i^*$可作为变量$v_i$的初始值用于精细化搜索（具体见模型求解部分）。最终得到的模型如下所示：
    \begin{align*}
            \max\limits_{x_i} T(x_i) = \sum_n O(t_n; x_i) 
            , \Delta t = 0.1s \\
            s.t.  \begin{cases}
               & v_i \in [70, 140], \\
               & t_{di} \leq t_{ei}, i=1,2,3. \\
               & \theta_i = \theta_{ik}^*, k=1,2,3,4,5. \\
            \end{cases}
    \end{align*}

    \subsubsection{模型求解}
    由问题二已求得无人机$FY_1$的最优投放策略，本问仅针对无人机$FY_2$和$FY_3$进行求解。代入初始数据
    $$
    u_2(0) = [12000, 1400, 1400], \quad u_3(0) = [6000, -3000, 700]
    $$
    根据航向选择模型，编写代码（见附件B）求解得到$FY_2$和$FY_3$的推荐航向和飞行速度如\cref{tab:pro4_dirc}：
    \begin{table}[!htbp]
    \caption{推荐航向和飞行速度}\label{tab:pro4_dirc} \centering
    \begin{tabular}{ccccc}
        \toprule[1.5pt]
        k &  $\theta^*_{2k}$ & $v^*_{2k}$ &  $\theta^*_{3k}$ & $v^*_{3k}$\\
        \midrule[1pt]
        1 & 322.30699062 & 133.417452 & 72.46982083 & 99.02928284\\
        2 & 322.24355395 & 126.175507 & 72.23633096 & 100.6299002\\
        3 & 322.14362733 & 122.8507499 & 72.1391272& 94.50079581\\
        4 & 322.17759766 & 111.1905521 & 72.19921093 & 107.5505793\\
        5 & 322.30699062 & 116.7179823 & 72.30440306 & 91.68939421\\
        \bottomrule[1.5pt]
    \end{tabular}
    \end{table}
 
    根据\cref{tab:pro4_dirc}，航向选择模型推荐的最佳航向均分布在一个较小的范围内，而飞行速度则差异较大。为保证结果的精确性，我们固定航向为$\theta^*_{21}$和$\theta^*_{31}$，并以推荐飞行速度$v^*_{21}$和$v^*_{31}$为初始值，对飞行速度、投放时间和起爆时间进行联合优化。
    优化算法与问题三一致，采用模拟退火算法（代码见附录B）。求解得到的投放策略结果如\cref{tab:pro4_result}所示:
    
    %\newpage

    \begin{table}[!htbp]
    \caption{问题四求解结果}\label{tab:pro4_result} \centering
    \begin{tabular}{c p{1.75cm} ccccc}
        \toprule[1.5pt]
        编号 &  航向角$\theta$ & 速度$v$ &  投放点坐标 & 起爆点坐标& 有效干扰时长& 干扰时间区间\\
         & $(^\circ)$ & $(m/s)$ &  &  & (s) &  \\
        \midrule[1pt]
        $FY_1$ & 352.82 & 99.99 & (17800, 0, 1800) & (17858, 7, 1798) & 4.60 & [2.55, 7.15]\\
        $FY_2$ & 322.31 & 129.98 & (13632, 139, 1400) & (13799, 10, 1387) & 3.28 & [17.50, 20.78]\\
        $FY_3$ & 72.47 & 100.01 & (6950, 9, 700)& (6966, 57, 699) &  3.56 & [32.06, 35.62]\\
        \bottomrule[1.5pt]
    \end{tabular}
\end{table}

    由\cref{tab:pro4_dirc}可知，三架无人机的投放策略均能有效遮蔽导弹$M_1$，且遮蔽时间区间均未重合，总有效遮蔽时长$T(x)$=11.44s，三架无人机投放的烟幕弹与导弹的局部运动轨迹如\cref{fig:pro4_3D}所示。
    \begin{figure}[!h]
    \centering
    \includegraphics[width=1.0\textwidth]{pro4_3D.png}
    \caption{问题四轨迹局部放大图}
    \label{fig:pro4_3D}
    \end{figure}

    \subsection{问题五模型的建立与求解}
    问题五为多架无人机协同对多个目标导弹进行干扰的场景进行建模与求解，由于决策变量较多，不妨限定一架无人机只针对一枚导弹进行遮蔽，将问题转化为多个“单无人机-单导弹”子问题进行处理。

    \subsubsection{模型建立}
    设无人机集合为 $\mathcal{F} = \{FY_1, FY_2, FY_3, FY_4, FY_5\}$，导弹集合为 $\mathcal{M} = \{M_1, M_2, M_3\}$。定义分配矩阵 $A_{i,j} \in \{0,1\}, 1 \leq i \leq 5, 1 \leq j \leq 3$，表示无人机 $FY_i$ 是否分配给导弹 $M_j$。
    每架无人机 $FY_i$ 的决策变量$x_i$为： 
    \begin{align*}
        & x_i = (\theta_i, v_i, td_{ik}, te_{ik}) \; \; 1 \leq i \leq 5, 1 \leq k \leq 3.
    \end{align*}
    其中$te_{ik}, td_{ik}$分别为第i架无人机投放的第k枚干扰弹的起爆时间和投放时间。

    对每颗导弹$M_j$，在任意时刻 $t$，其遮蔽状态函数 $O_M(t; x)$ 定义为：
    \begin{align*}
        O_M(t; x) =
        \begin{cases}
            1, & \exists i \text{ s.t. } A_{i,j} = 1 \text{ and } \max d_{PM \rightarrow C}(t) \leq 10, \\
            0, & \text{otherwise}.
        \end{cases}
    \end{align*}

    在多枚干扰弹遮蔽取并集情况下，导弹 $M$ 的总遮蔽时间 $T_M$ 计算为：
    \begin{align*}
        T_M(x) = \sum_{t=0}^{t_{\text{flight}}}O_M(t; x) \Delta t, 
    \end{align*}
    其中$t_{\text{flight}}$为最长导弹飞行总时长，根据题干条件为导弹$M_1$飞行时长66.67s。

    为尽可能均衡分配无人机对不同导弹进行遮蔽，我们取各导弹总遮蔽时长中最小值，将该值的最大化作为\textbf{目标函数:}
    \[
    \max \min_{M \in \mathcal{M}} T_M(x)
    \]

    根据题干所给约束条件，问题五模型建立如下：
    \begin{align*}
        \max \min_{M \in \mathcal{M}} T_M(x) = \max \min_{M \in \mathcal{M}} \sum_{t=0}^{t_{\text{flight}}}O_M(t; x) \Delta t, \\
        s.t. 
        \begin{cases}
        & A_{i,j} \in \{0,1\}, \\
         & v_i \in [70, 140], \\
        & \theta_i \in [0, 2\pi), \\
        & td_{i,k} \leq te_{i,k}, \; k=1,2,3, \\
        & td_{i,k+1} - td_{i,k} \geq 1, \; k=1,2,3 \\
        \end{cases}
    \end{align*}

\subsubsection{模型求解}

\textbf{第一层：单无人机-单导弹参数优化}

    对每一架无人机$FY_i$，我们分别令其针对$M_j, j=1,2,3$进行遮蔽，即分别令$A_{i,j}=1, j=1,2,3$.给定任意$(FY_i, M_j)$组合，使用模拟退火算法结合网格搜索，为其寻找最优投放策略$x_i$，使得$T_{M_j}$最大。

    我们借鉴问题四的思路，对飞行角度$\theta$、飞行速度$v$、投放时间$td$和起爆时间$te$进行分层优化。第一轮搜索中，我们固定飞行速度$v$，对飞行角度$\theta$进行搜索，使用模拟退火算法求解每组固定$(\theta, v)$条件下最佳遮蔽时间以及$td$、$te$投放策略，
    对每个$(FY_i, M_j)$给出其最优飞行角度如\cref{tab:最优角度}。

    在给定最优飞行角度后，我们对$v$进行变步长网格搜索，同样使用模拟退火算法对投放时间、起爆时间进行优化，最终确定每个$(FY_i, M_j)$组合的最优投放策略。
    \begin{table}[htbp]
    \centering
    \caption{角度优化结果}
    \label{tab:最优角度}
    \begin{tabular}{|c|c|c|c|c|c|}
    \hline
    \diagbox{导弹}{无人机} & $FY_1$& $FY_2$ & $FY_3$ & $FY_4$ & $FY_5$ \\
    \hline
    $M_1$ & 180 & 239 & 78 & 350 & 106\\
    \hline
    $M_2$ & 350 & 286 & 88 & 350 & 119\\
    \hline
    $M_3$ & 350 & 231 & 116 & 350 & 121\\
    \hline
    \end{tabular}
    \end{table}

    

\textbf{第二层：无人机-导弹分配方案优化}

    通过第一层的优化，我们得到每一对$(FY_i, M_j)$对应的最优$T_{M_j}$.理论上$(FY_i, M_j)$所组成的分配矩阵$A$共有 $3^5 = 243$种可能。我们遍历每种分配方案，并通过以下流程进行评估：
\begin{enumerate}
    \item 收集爆炸时间：获取分配给每枚导弹的所有干扰弹爆炸时间
    \item 计算并集遮蔽时间：使用时间离散化方法，计算每枚导弹的遮蔽时间并集
    \item 评估方案质量：以最小遮蔽时间为目标函数值
\end{enumerate}



\subsubsection{求解结果与分析}

\textbf{最优分配方案}
遍历所有分配方案，得到最优分配结果及其对应飞行策略如\cref{tab:pro5_result}所示：

\begin{table}[h!]
    \centering
    \caption{问题五投放策略结果}
    \label{tab:pro5_result}
    \begin{tabular}{lcccc}
        \toprule
        无人机 & 分配导弹 & 飞行角度$\theta$(°) & 飞行速度$v$(m/s) & 遮蔽时间(s) \\
        \midrule
        FY1 & M1 & 180 & 120.8 & 6.04 \\
        FY2 & M2 & 286 & 133.0 & 10.98 \\
        FY3 & M3 & 116 & 140.0 & 2.63 \\
        FY4 & M1 & 350 & 120.3 & 0.00 \\
        FY5 & M3 & 106 & 138.7 & 5.08 \\
        \bottomrule
    \end{tabular}
\end{table}
此时三枚导弹总体被遮蔽时间$\sum\limits_{M \in \mathcal{M}}T_M = 24.73s$，各个烟幕弹的投放坐标、起爆坐标、遮蔽导弹和遮蔽时长如\cref{tab:pro5_烟幕弹}所示，
各个烟幕弹的投放时间、起爆时间和导弹被有效遮蔽的时间区间如\cref{fig:pro5_result}所示：
\begin{table}[h!]
    \centering
    \caption{问题五各烟幕弹数据}
    \label{tab:pro5_烟幕弹}
    \begin{tabular}{lccccc}
        \toprule
        无人机& 烟幕弹编号 & 导弹 & 投放坐标 & 起爆坐标 & 遮蔽时间(s) \\
        \midrule
         & 1& M1 & (17800,0,1800) & (17403,0, 1747) & 2.95 \\
        FY1 & 2& M1 & (173600,1800) & (16753, 0, 1676) & 2.24 \\
         & 3& M1 & (17146, 0,1800)  & (16486, 0, 1654) & 0.93 \\
         \hline
         & 1& M2 & (12179,775,1400) & (12277, 432, 1365) & 3.84 \\
        FY2 & 2& M2 & (12216,647,1400) & (12283, 413, 1384) & 3.77 \\
         & 3& M2 & (12268, 464, 1400) & (12286, 400, 1399) & 3.37 \\
         \hline
         & 1& M3 & (4958, -863, 700) & (4558, -44, 492) & 2.63 \\
        FY3 & 2& M3 & (4589, -108, 700) & (4375, 332, 640) & 0 \\
         & 3& M3 & (4221, 647, 700) & (4006, 1087, 640) & 0 \\
         \hline
         & 1& M1 & (11177, 1968, 1800) & (12362, 1759, 1310) & 0 \\
        FY4 & 2& M1 &(12599, 1717, 1800) & (13025, 1642, 1736)& 0 \\
         & 3& M1 & (11059, 1989, 1800)& (11592, 1895, 1700) & 0 \\
         \hline
         & 1& M3 & (12192, -656, 1300) & (12028, -382, 1274) & 3.61 \\
        FY5 & 2& M3 & (12117,-531, 1300) & (12039, -400, 1294) & 1.47 \\
         & 3& M3 & (11689, 182, 1300) & (11439, 598, 1240) & 0 \\
        \bottomrule
    \end{tabular}
\end{table}


\begin{figure}[!h]
    \centering
    \includegraphics[width=1.0\textwidth]{pro5.png}
    \caption{投放/起爆时间和遮蔽区间示意图}
    \label{fig:pro5_result}
\end{figure}

\cref{fig:pro5_result}中每架无人机横线上烟幕弹投放时间和起爆时间分别用圆圈和叉型标记表示，同一颜色的两点对应同一枚烟幕弹；导弹横线涂蓝色部分表示其被有效遮蔽时间区间。

\section{模型检验与评价}

\subsection{模型检验}
    \subsubsection{问题二模型检验}
    由于PSO优化算法中可能存在局部最优解导致的误差，我们经如下两步对问题二的模型进行检验：
    \begin{enumerate}
        \item 随机生成可行解并评估:我们在问题二给定的题干约束范围下随机生成了100组可行解，并调整最优解附近生成的可行解数目更加密集。对每个可行解分别计算有效遮蔽时间，
    与最优解进行比较，结果如\cref{fig:可行解对比}所示。
        \item 粒子群算法灵敏度分析:在粒子群优化(PSO)算法中，其求解质量高度依赖认知系数和社会系数的设置，
    我们取认知系数$c_1 \in [1.5,2.0]$，社会系数$c_2 \in [2.0,2.5]$，并以步长0.1取值，代入问题二模型分别求解，观察不同参数组合下
    遮蔽时间的情况。最终结果如3D曲面图\cref{fig:curve}所示。
    \end{enumerate}

    \begin{figure}
    \centering
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=.8\textwidth]{check.png}
        \subcaption{随机可行解对比图}
        \label{fig:可行解对比}
    \end{minipage}
    \begin{minipage}[c]{0.48\textwidth}
        \centering
        \includegraphics[width=.8\textwidth]{pso_sensitivity.png}
        \subcaption{灵敏度分析曲面图}
        \label{fig:curve}
    \end{minipage}
    \caption{问题二模型检验示意图}
    \label{fig:pro2_check}
    \end{figure}

    观察图像，我们发现\cref{fig:可行解对比}中所有可行解均不高于问题二模型优化所得到的最优解；由\cref{fig:curve}可知认知系数和社会系数对PSO算法的求解结果有较大影响，
    不同参数组合所得到的最优解大致分布在[3.4, 4.6]间，且图像在$c_1 = 2.0, c_2 = 2.0$点处出现明显高峰，证明我们在问题二中的参数设置为最优。
    
    \subsubsection{问题四模型检验}
    问题四相比前序问题的核心变化在于航迹选择模型的使用。$FY_2$和$FY_3$根据航迹选择模型的推荐航向和飞行速度进行优化，而
    $FY_1$的航向和速度沿用问题二所求得的最优策略。由于问题二可视为问题四的子问题，因此本部分借用问题二情境验证航迹选择模型的准确性。

    我们将问题二中$FY_1$和$M_1$的参数根据代入航迹选择模型，通过模型中的评估函数得到评分最高的3组飞行角度及其对应预估速度，并与问题二中的结果进行对比。如\cref{tab:pro4_check}所示。
    \begin{table}[htbp]
    \centering
    \caption{航迹选择模型检验}
    \label{tab:pro4_check}
    \begin{tabular}{|c|c|c|c|c|}
    \hline
     \diagbox{求解结果}{模型} & \multicolumn{1}{c}{} & \multicolumn{1}{c}{航迹选择模型} &  & 问题二模型 \\
        \hline
                编号          & 1& 2 & 3 & $FY_1$ \\
    \hline
    飞行角度$\theta(^\circ)$ & 359.59 & 359.55 & 359.63 & 352.82 \\
    \hline
    速度$v(m/s)$ & 98.45 & 78.87 & 131.11 & 99.99\\
    \hline
    \end{tabular}
    \end{table}


    \cref{tab:pro4_check}中航迹选择模型所得到的飞行方向和速度如前三列所示：$FY_1$的数据为问题二模型优化得到的，。由表中结果与$FY_1$用PSO算法得出的数据，我们发现，航迹选择模型对$FY_1$的优化结果与PSO算法
    的优化结果略有差异，但是差异不大，特别是第一组数据，在速度上非常接近，飞行角度虽略有差异，但是
    也在可接受范围之内。因此，我们认为本题中使用的航迹选择模型是合理的。

    \subsection{模型评价}
        \subsubsection{模型优点}
        \begin{enumerate}
            \item \textbf{模型间关联性强：}问题三模型经过适当推广和参数修改可作为问题四的子模型，同时经过验证，问题四的模型在问题二情境下求解结果与问题二自身模型保持高度的一致性，说明我们的各个模型之间思路连贯，联系紧密。
            \item \textbf{模型求解方法多样：}我们针对不同问题的特点，选择了粒子群、模拟退火算法、变步长网格搜索等多种方式进行求解，并通过参数调节、结果对比等方式保证了求解结果的可靠性。
        \end{enumerate}
        \subsubsection{模型缺点}
        \begin{enumerate}
            \item \textbf{模型的解有概率并非全局最优：}在模型求解过程中，我们使用粒子群算法、模拟退火算法等启发式算法，可能由于初始参数的设置陷入局部最优。我们通过调整模型参数，细化搜索步长等方式尽可能避免此问题。
            \item \textbf{问题五未考虑一架无人机对多个导弹进行遮蔽：} 最初我们考虑该情况，并采用遗传算法进行求解，但由于模型过于复杂，求解时间代价高，我们基于问题四的分析，将问题五简化为多个“单无人机-单导弹”子问题进行处理，可能会损失部分最优解。
            \item \textbf{未考虑多个烟幕弹共同遮蔽目标的情况：}在实际场景中，单独的一个烟幕云团可能无法对目标实现完整遮蔽，需要多个云团的共同配合。在本模型建立过程中为简化求解运算忽略了这一情形，可能会影响最终的遮蔽效果。
        \end{enumerate}

%参考文献
\begin{thebibliography}{9}%宽度9
    \bibitem[1]{ref1}
    罗瑞耀,王得霖,罗威,等.烟幕弹应对察打一体无人机的投放策略研究\allowbreak[J]..光电技术应用,2022,37(06):90-98.
    \bibitem[2]{ref2}商家鹏.基于改进粒子群算法的物联网链路负载均衡优化方法[J].电脑编程技巧与维护,2025,(08):28-30+52.DOI:10.16184/j.cnki.comprg.2025.08.001.
    \bibitem[3]{ref3}潘翔,徐凯.基于遗传模拟退火算法的智能仓储多AGV调度研究[J/OL].浙江工业大学学报,1-7[2025-09-07].https://link.cnki.net/urlid/33.1193.T.20250808.0941.002.
    \bibitem[4]{ref4}陆青梅,赵山林,高媛.基于模拟退火遗传算法的舰船编队网络优化调度方法[J].舰船科学技术,2025,47(10):155-160.
    \bibitem[5]{ref5}时光志.基于遗传模拟退火算法的多船避碰决策[J].舰船科学技术,2025,47(09):154-159.
    \bibitem[6]{ref6}王翰章,张雪涛,刘一莎,等.基于子母式多无人机协同的多目标追踪[J].机器人,2025,47(03):348-360.DOI:10.13973/j.cnki.robot.240352.
    \bibitem[7]{ref7}Hu Y ,Zhang X ,Deng K , et al.Simulated annealing-based nonmonotone conjugate gradient method for unconstrained optimization with applications[J].Optimization and Engineering,2025,(prepublish):1-42.
    \bibitem[8]{ref8}Ma B ,Ji Y ,Fang L .A Multi-UAV Formation Obstacle Avoidance Method Combined with Improved Simulated Annealing and an Adaptive Artificial Potential Field[J].Drones,2025,9(6):390-390.
    \bibitem[9]{ref9}Zhuo M ,Feng Y ,Yang P , et al.Optimizing Topology in Satellite–UAV Collaborative IoT: A Graph Partitioning Simulated Annealing Approach[J].Drones,2024,8(2):44-.
    \bibitem[10]{ref10}Cong C ,Yibai L ,Guangqiao C , et al.Research on Dynamic Scheduling Model of Plant Protection UAV Based on Levy Simulated Annealing Algorithm[J].Sustainability,2023,15(3):1772-1772.
    \bibitem[11]{ref11}Wang S ,Qiao P ,Yue Q , et al.Research on Dynamic Particle Swarm Optimization for Multi-Objective Reconnaissance Task Allocation of UAVs[J].Drones,2025,9(8):556-556.
\end{thebibliography}


\newpage
%附录
\begin{appendices}
 \section{支撑材料文件列表}
    \begin{table}[H]
        \centering
        \caption{支撑材料文件列表}
        \begin{tabular}{clc}
            \toprule[1.5pt]
            序号 & 文件名 & 说明 \\
            \midrule[1pt]
            1 & problem1.py & 问题一求解代码 \\
            
            \bottomrule[1.5pt]
        \end{tabular}
    \end{table}
\section{代码清单}
\noindent\textbf{1.问题一求解代码(problem1.py)}
\begin{lstlisting}[language=python, numbers=left]
import math
import xlsxwriter

# 存储xlsx文件
workbook = xlsxwriter.Workbook('problem1.xlsx')
worksheet = workbook.add_worksheet()
worksheet.write(0, 0, 'Time (s)')
worksheet.write(0, 1, 'Max_d2')

workbook2 = xlsxwriter.Workbook('problem1_pos.xlsx')
worksheet2 = workbook2.add_worksheet()
worksheet2.write(0, 1, 'C_x')
worksheet2.write(0, 2, 'C_y')
worksheet2.write(0, 3, 'C_z')
worksheet2.write(0, 4, 'M_x')
worksheet2.write(0, 5, 'M_y')
worksheet2.write(0, 6, 'M_z')

FY1 = (17800.0, 0.0, 1800.0)    # 无人机 FY1 初始位置
T = (0.0, 200.0, 0.0)           # 真目标中心
M1 = (20000.0, 0.0, 2000.0)     # 导弹 M1 初始位置
v_uav = 120.0                   # UAV 速度
t_drop = 1.5                    # 投放时间
fuse = 3.6                      # 引信延时
g = 9.8                         # 重力加速度
cloud_r = 10.0                  # 云团半径
cloud_sink = 3.0                # 云团下沉速度
v_missile = 300.0               # 导弹速度

# FY1飞向原点的水平单位向量
hx, hy = -FY1[0], -FY1[1]
h_norm = math.hypot(hx, hy)
ux, uy = hx / h_norm, hy / h_norm

# 起爆时刻与位置
t_explosion = t_drop + fuse
C0 = (
    FY1[0] + v_uav * ux * t_explosion,   # x
    FY1[1] + v_uav * uy * t_explosion,   # y
    FY1[2] - 0.5 * g * fuse**2           # z
)

# 导弹速度向量
to_origin = (-M1[0], -M1[1], -M1[2])
norm = math.sqrt(sum(c*c for c in to_origin))
v_m = tuple(c / norm * v_missile for c in to_origin)

# 真目标圆柱采样
r = 7
h = 10
N_THETA = 80            # 周向采样数（可调：16~48）
N_Z = 300                 # 轴向采样层数（可调：3~7）

xc, yc, zc = T
theta_vals = [2*math.pi*k/N_THETA for k in range(N_THETA)]
z_vals = [zc + k*(h/(N_Z-1)) for k in range(N_Z)]

# 在圆柱侧面（各层圆周）生成采样点；可选再加轴线点提升覆盖
target_samples = []
for z in z_vals:
    for th in theta_vals:
        sample = (
            xc + r * math.cos(th),
            yc + r * math.sin(th),
            z
        )
        target_samples.append(sample)


# 采样计算遮蔽区间
dt = 0.01
shield_flags = []
times = [i*dt for i in range(int(20/dt)+1)]

max_d2_list = []
pos_list = []
for tau in times:
    # 导弹位置
    t_global = t_explosion + tau
    Mx = M1[0] + v_m[0]*t_global
    My = M1[1] + v_m[1]*t_global
    Mz = M1[2] + v_m[2]*t_global

    # 云团中心
    Cx, Cy, Cz = C0[0], C0[1], C0[2] - cloud_sink*tau
    if tau % 0.5 == 0:
        print("云团：", Cx, Cy, Cz, "导弹：", Mx, My, Mz)
        pos_list.append((tau, Cx, Cy, Cz, Mx, My, Mz))
    # 判断云团球体是否与 M→T 线段相交
    max_d2 = -float('inf')
    for Px, Py, Pz in target_samples:
        # 计算点到线段距离平方
        vx, vy, vz = Px - Mx, Py - My, Pz - Mz   # M->P
        wx, wy, wz = Cx - Mx, Cy - My, Cz - Mz   # M->C
        seg_len2 = vx*vx + vy*vy + vz*vz

        if seg_len2 == 0.0:
            d2 = wx*wx + wy*wy + wz*wz
        else:
            u = (wx*vx + wy*vy + wz*vz) / seg_len2
            if u <= 0.0:
                # 最近点是 M
                dx, dy, dz = wx, wy, wz
            elif u >= 1.0:
                # 最近点是 P
                dx, dy, dz = Cx - Px, Cy - Py, Cz - Pz
            else:
                # 最近点在 M-P 之间
                projx = Mx + u*vx
                projy = My + u*vy
                projz = Mz + u*vz
                dx, dy, dz = Cx - projx, Cy - projy, Cz - projz
            d2 = dx*dx + dy*dy + dz*dz

        if d2 > max_d2:
            max_d2 = d2

    max_d2_list.append(math.sqrt(max_d2))
    shield_flags.append(max_d2 <= cloud_r ** 2)


# 提取连续遮蔽区间
intervals, total_time = [], 0.0
i = 0
while i < len(times):
    if shield_flags[i]:
        j = i
        while j+1 < len(times) and shield_flags[j+1]:
            j += 1
        start, end = times[i], times[j]+dt
        intervals.append((start, end))
        total_time += end - start
        i = j + 1
    else:
        i += 1

print("Explosion time:", t_explosion, "s")
print("Explosion position C0:", C0)
print("有效遮蔽区间 (s after explosion):", intervals)
print("总有效遮蔽时长:", total_time, "s")

for i, (time, max_d2) in enumerate(zip(times, max_d2_list)):
    worksheet.write(i + 1, 0, time)
    worksheet.write(i + 1, 1, max_d2)

for i, (tau, Cx, Cy, Cz, Mx, My, Mz) in enumerate(pos_list):
    worksheet2.write(i + 1, 0, tau)
    worksheet2.write(i + 1, 1, Cx)
    worksheet2.write(i + 1, 2, Cy)
    worksheet2.write(i + 1, 3, Cz)
    worksheet2.write(i + 1, 4, Mx)
    worksheet2.write(i + 1, 5, My)
    worksheet2.write(i + 1, 6, Mz)

workbook.close()
workbook2.close()
print("Data saved to 'problem1.xlsx'")
print("Data saved to 'problem1_pos.xlsx'")
\end{lstlisting}
\end{appendices}

\end{document} 